#####################
# Sleep Project - Pedro Bessone, Gautam Rao, Heather Schofield, Frank Schilbach, and Mattie Toma
# Purpose: Provides adjusted p-values for table 17 from the Appendix (Main Treatment Effects, Pooling Night-Sleep Treatments and Including Multiple Hypothesis Testing Corrections). 
#          Takes the unadjusted p-values from table_A17.do and performs the
#          Hochberge corrections, and exports the adjusted p-values.
# Last edited: 07 May 2021
#####################

# install.packages("BiocManager")
require(multtest)
require(rio)
require(tidyverse)

rm(list = ls())

options(scipen = 99)

df.p = import("Datasets/pvals_unadjusted_table_instrumental.csv") %>% 
  as_tibble()
colnames(df.p) = colnames(df.p) %>% str_remove('[0-9]')

df.family = df.p %>% 
  filter( var_name %in%  c("earnings", "wellbeing", "cogindex", "pref") )

df.work = df.p %>% 
  filter( var_name %in%  c("prod", "typing", "output"))

df.wb = df.p %>% 
  filter( var_name %in%  c("physical", "psych"))

df.cog = df.p %>% 
  filter( var_name %in%  c("cogfunction", "attention"))

df.pref = df.p %>% 
  filter( var_name %in%  c("timepref", "social", "risk"))

# Parameters
# corrections
corrections = c("BH", "BY", "hochberg", "holm", "TSBH", "hommel")
alpha = 0.05 #(for rejection threshold when using the BKY06 method)


p.adjust.vec = function(p.unadjust, corrections, alpha){
  p.adj.mat = matrix(NA, 
                     ncol = length(corrections), 
                     nrow = length(p.unadjust))
  colnames(p.adj.mat) = corrections
  
  for(i in seq(1, length(corrections))){
    corr = corrections[i]
    
    if(corr == "TSBH"){
      ls.p.adj = multtest::mt.rawp2adjp(p.unadjust,
                                        proc = corr,
                                        alpha = 0.05)  
      p.adj.mat[,i] = ls.p.adj$adjp[,2]
    } else{
      ls.p.adj = p.adjust(p.unadjust, method = corr )
      p.adj.mat[,i] = ls.p.adj
    }
    
    
  }
  
  return(p.adj.mat)
} 

# Comparison: NS vs. Control

df.family.ns.cont = df.family %>% 
  select(var_name, p_val_ns) %>% 
  mutate(comparison = "NS vs. Control") %>% 
  rename(p_val_unadjust = p_val_ns) %>% 
  arrange(p_val_unadjust) 

p.unadjust = df.family.ns.cont$p_val_unadjust

p.adj.mat = p.adjust.vec(p.unadjust, corrections = corrections, alpha = alpha) %>% as_tibble()

df.family.ns.cont = df.family.ns.cont %>% 
  bind_cols(p.adj.mat)

# Comparison: Nap vs. Control

df.family.nap.cont = df.family %>% 
  select(var_name, p_val_nap) %>% 
  mutate(comparison = "Nap vs. Control") %>% 
  rename(p_val_unadjust = p_val_nap) %>% 
  arrange(p_val_unadjust)

p.unadjust = df.family.nap.cont$p_val_unadjust

p.adj.mat = p.adjust.vec(p.unadjust, corrections = corrections, alpha = alpha) %>% as_tibble()

df.family.nap.cont = df.family.nap.cont %>% 
  bind_cols(p.adj.mat)

# Comparison: Nap vs. NS

df.family.nap.ns = df.family %>% 
  select(var_name, p_val_d_nap_ns) %>% 
  mutate(comparison = "Nap vs. NS") %>% 
  rename(p_val_unadjust = p_val_d_nap_ns) %>% 
  arrange(p_val_unadjust)

p.unadjust = df.family.nap.ns$p_val_unadjust

p.adj.mat = p.adjust.vec(p.unadjust, corrections = corrections, alpha = alpha) %>% as_tibble()

df.family.nap.ns = df.family.nap.ns %>% 
  bind_cols(p.adj.mat)

df.family.all = df.family.ns.cont %>% 
  bind_rows(df.family.nap.cont) %>% 
  bind_rows(df.family.nap.ns) %>% 
  relocate(var_name, comparison)

df.family.all = df.family.all %>% 
  rename(BKY06 = TSBH) %>% 
  arrange(comparison, var_name)

# Work Family


# Comparison: NS vs. Control

df.work.ns.cont = df.work %>% 
  select(var_name, p_val_ns) %>% 
  mutate(comparison = "NS vs. Control") %>% 
  rename(p_val_unadjust = p_val_ns) %>% 
  arrange(p_val_unadjust)

p.unadjust = df.work.ns.cont$p_val_unadjust

p.adj.mat = p.adjust.vec(p.unadjust, corrections = corrections, alpha = alpha) %>% as_tibble()

df.work.ns.cont = df.work.ns.cont %>% 
  bind_cols(p.adj.mat)

# Comparison: Nap vs. Control

df.work.nap.cont = df.work %>% 
  select(var_name, p_val_nap) %>% 
  mutate(comparison = "Nap vs. Control") %>% 
  rename(p_val_unadjust = p_val_nap) %>% 
  arrange(p_val_unadjust)

p.unadjust = df.work.nap.cont$p_val_unadjust

p.adj.mat = p.adjust.vec(p.unadjust, corrections = corrections, alpha = alpha) %>% as_tibble()

df.work.nap.cont = df.work.nap.cont %>% 
  bind_cols(p.adj.mat)

# Comparison: Nap vs. NS

df.work.nap.ns = df.work %>% 
  select(var_name, p_val_d_nap_ns) %>% 
  mutate(comparison = "Nap vs. NS") %>% 
  rename(p_val_unadjust = p_val_d_nap_ns)  %>% 
  arrange(p_val_unadjust)

p.unadjust = df.work.nap.ns$p_val_unadjust

p.adj.mat = p.adjust.vec(p.unadjust, corrections = corrections, alpha = alpha) %>% as_tibble()

df.work.nap.ns = df.work.nap.ns %>% 
  bind_cols(p.adj.mat)


df.work.all = df.work.ns.cont %>% 
  bind_rows(df.work.nap.cont) %>% 
  bind_rows(df.work.nap.ns) %>% 
  relocate(var_name, comparison)

df.work.all = df.work.all %>% 
  rename(BKY06 = TSBH) %>% 
  arrange(comparison, var_name)

# Well-Being Family


# Comparison: NS vs. Control

df.wb.ns.cont = df.wb %>% 
  select(var_name, p_val_ns) %>% 
  mutate(comparison = "NS vs. Control") %>% 
  rename(p_val_unadjust = p_val_ns) %>% 
  arrange(p_val_unadjust)

p.unadjust = df.wb.ns.cont$p_val_unadjust

p.adj.mat = p.adjust.vec(p.unadjust, corrections = corrections, alpha = alpha) %>% as_tibble()

df.wb.ns.cont = df.wb.ns.cont %>% 
  bind_cols(p.adj.mat)

# Comparison: Nap vs. Control

df.wb.nap.cont = df.wb %>% 
  select(var_name, p_val_nap) %>% 
  mutate(comparison = "Nap vs. Control") %>% 
  rename(p_val_unadjust = p_val_nap) %>% 
  arrange(p_val_unadjust)

p.unadjust = df.wb.nap.cont$p_val_unadjust

p.adj.mat = p.adjust.vec(p.unadjust, corrections = corrections, alpha = alpha) %>% as_tibble()

df.wb.nap.cont = df.wb.nap.cont %>% 
  bind_cols(p.adj.mat)

# Comparison: Nap vs. NS

df.wb.nap.ns = df.wb %>% 
  select(var_name, p_val_d_nap_ns) %>% 
  mutate(comparison = "Nap vs. NS") %>% 
  rename(p_val_unadjust = p_val_d_nap_ns)  %>% 
  arrange(p_val_unadjust)

p.unadjust = df.wb.nap.ns$p_val_unadjust

p.adj.mat = p.adjust.vec(p.unadjust, corrections = corrections, alpha = alpha) %>% as_tibble()

df.wb.nap.ns = df.wb.nap.ns %>% 
  bind_cols(p.adj.mat)


df.wb.all = df.wb.ns.cont %>% 
  bind_rows(df.wb.nap.cont) %>% 
  bind_rows(df.wb.nap.ns) %>% 
  relocate(var_name, comparison)

df.wb.all = df.wb.all %>% 
  rename(BKY06 = TSBH) %>% 
  arrange(comparison, var_name)

# Cognition Family


# Comparison: NS vs. Control

df.cog.ns.cont = df.cog %>% 
  select(var_name, p_val_ns) %>% 
  mutate(comparison = "NS vs. Control") %>% 
  rename(p_val_unadjust = p_val_ns) %>% 
  arrange(p_val_unadjust)

p.unadjust = df.cog.ns.cont$p_val_unadjust

p.adj.mat = p.adjust.vec(p.unadjust, corrections = corrections, alpha = alpha) %>% as_tibble()

df.cog.ns.cont = df.cog.ns.cont %>% 
  bind_cols(p.adj.mat)

# Comparison: Nap vs. Control

df.cog.nap.cont = df.cog %>% 
  select(var_name, p_val_nap) %>% 
  mutate(comparison = "Nap vs. Control") %>% 
  rename(p_val_unadjust = p_val_nap) %>% 
  arrange(p_val_unadjust)

p.unadjust = df.cog.nap.cont$p_val_unadjust

p.adj.mat = p.adjust.vec(p.unadjust, corrections = corrections, alpha = alpha) %>% as_tibble()

df.cog.nap.cont = df.cog.nap.cont %>% 
  bind_cols(p.adj.mat)

# Comparison: Nap vs. NS

df.cog.nap.ns = df.cog %>% 
  select(var_name, p_val_d_nap_ns) %>% 
  mutate(comparison = "Nap vs. NS") %>% 
  rename(p_val_unadjust = p_val_d_nap_ns) %>% 
  arrange(p_val_unadjust)

p.unadjust = df.cog.nap.ns$p_val_unadjust

p.adj.mat = p.adjust.vec(p.unadjust, corrections = corrections, alpha = alpha) %>% as_tibble()

df.cog.nap.ns = df.cog.nap.ns %>% 
  bind_cols(p.adj.mat)


df.cog.all = df.cog.ns.cont %>% 
  bind_rows(df.cog.nap.cont) %>% 
  bind_rows(df.cog.nap.ns) %>% 
  relocate(var_name, comparison)

df.cog.all = df.cog.all %>% 
  rename(BKY06 = TSBH) %>% 
  arrange(comparison, var_name)

# Preferences Family


# Comparison: NS vs. Control

df.pref.ns.cont = df.pref %>% 
  select(var_name, p_val_ns) %>% 
  mutate(comparison = "NS vs. Control") %>% 
  rename(p_val_unadjust = p_val_ns) %>% 
  arrange(p_val_unadjust)

p.unadjust = df.pref.ns.cont$p_val_unadjust

p.adj.mat = p.adjust.vec(p.unadjust, corrections = corrections, alpha = alpha) %>% as_tibble()

df.pref.ns.cont = df.pref.ns.cont %>% 
  bind_cols(p.adj.mat)

# Comparison: Nap vs. Control

df.pref.nap.cont = df.pref %>% 
  select(var_name, p_val_nap) %>% 
  mutate(comparison = "Nap vs. Control") %>% 
  rename(p_val_unadjust = p_val_nap) %>% 
  arrange(p_val_unadjust)

p.unadjust = df.pref.nap.cont$p_val_unadjust

p.adj.mat = p.adjust.vec(p.unadjust, corrections = corrections, alpha = alpha) %>% as_tibble()

df.pref.nap.cont = df.pref.nap.cont %>% 
  bind_cols(p.adj.mat)

# Comparison: Nap vs. NS

df.pref.nap.ns = df.pref %>% 
  select(var_name, p_val_d_nap_ns) %>% 
  mutate(comparison = "Nap vs. NS") %>% 
  rename(p_val_unadjust = p_val_d_nap_ns) %>% 
  arrange(p_val_unadjust)

p.unadjust = df.pref.nap.ns$p_val_unadjust

p.adj.mat = p.adjust.vec(p.unadjust, corrections = corrections, alpha = alpha) %>% as_tibble()

df.pref.nap.ns = df.pref.nap.ns %>% 
  bind_cols(p.adj.mat)

df.pref.all = df.pref.ns.cont %>% 
  bind_rows(df.pref.nap.cont) %>% 
  bind_rows(df.pref.nap.ns) %>% 
  relocate(var_name, comparison)

df.pref.all = df.pref.all %>% 
  rename(BKY06 = TSBH) %>% 
  arrange(comparison, var_name)

# Combine all


df.all = bind_rows(
  df.family.all,
  df.work.all,
  df.wb.all,
  df.cog.all,
  df.pref.all
) %>% 
  mutate(
    family = rep(
      c("across-family", "work", "WB", "cog", "preferences"),
      c(4*3,3*3,2*3, 2*3, 3*3)
    )
  ) 

export(df.all, "Datasets/pval_adj_all_instrumental.xlsx")

